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ABSTRACT 


The  effects  of  combined  radiative  and  convective  heat  transfer  from 
arrays  of  longitudinal  fins  of  arbitrary  profile  are  analyzed  subject  to 
non-uniform  surface  emissivity  and  non-uniform  surface  film  coeffi¬ 
cients.  Consideration  is  given  to  radiative  interactions  between  adja¬ 
cent  fins  and  between  fins  and  the  base  surface.  Solution  of  the  defining 
differential  equation  for  fin  temperature  distribution  is  obtained  through 
an  iterative  application  of  the  B.  G.  Galerkin  variational  technique. 
Application  of  the  method  of  solution  is  made  to  fins  of  parabolic,  tri¬ 
angular,  and  inverse  parabolic  profile  subject  solely  to  the  radiative 
mode  of  heat  transfer.  Effects  in  variations  of  the  dimensionless 
radiation  number.  Nr,  fin  spacing,  S,  and  fin  surface  emissivity,  e, 
are  investigated.  Findings  of  the  study  reveal  that  for  the  pure  radia¬ 
tive  mode,  fins  can  enhance  the  heat  transfer  between  the  base  and  the 
surroundings  only  for  the  case  of  low  fin  and  base  surface  emissivity. 
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CHAPTER  I 

INTRODUCTION 

Heat  transfer  phenomena  associated  with  convection 
from  arrays  of  longitudinal . fins  has  been  the  subject  of 
many  theoretical  and  experimental  investigations  since  early 
in  this  century.  Most  early  investigators  (1,  2,  3)^ 
limited  their  interest  to  single  fins  of  simplified  geometry 
and  employed  certain  classical  assumptions  such  as  uniform 
convective  film  coefficients  and  negligible  influence  of  the 
fin  surface  slope.  In  all  cases  the  influence  of  radiation 
to  the  environment  and  to  surrounding  surfaces  was  neglected. 

In  1951,  Ghai  (4)  experimentally  demonstrated  that 
for  certain  fin  configurations  the  convective  film  coeffi¬ 
cients  were  larger  near  the  tip  than  near  the  fin  base. 

Ghai  concluded  that  the  classical  assumption  of  uniform  film 
coefficients  was  not,  in  general,  valid. 

Eraslan  and  Frost  (5)  have  solved  the  problem  of  the 
longitudinal  convection  fin  subject  to  non-uniform  film 
coefficients.  Their  study  indicates  that  fins  previously 
considered  as  being  optimum  when  analyzed  assuming  uniform 
film  coefficients,  did  not  necessarily  remain  so  when  the 

^Numbers  in  parentheses  refer  to  similarly  numbered 
references  in  the  bibliography. 
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effects  of  non-uniform  film  coefficients  were  considered. 
Eraslan  and  Frost  applied  a  variational  method  due  to  B.  G. 
Galerkin  in  their  solution. 

With  the  advent  of  the  space  age,  the  pure  radiation 
fin  became  of  significance.  Several  authors  (6,  7)  consid¬ 
ered  a  single  fin  radiating  to  an  isothermal  environment 
where  interactions  were  neglected.  Eckert,  Irvine,  and 
Sparrow  (8,  9),  in  1960,  published  a  series  of  papers  which 
formulated  the  basic  equations  defining  radiative  inter¬ 
action  in  arrays  of  longitudinal  fins.  Their  problem 
considered  radiative  interactions  of  the  fins  with  adjacent 
fins  and  with  base  surfaces.  These  authors,  however, 
limited  the  application  of  their  analysis  to  fairly  simpli¬ 
fied  shapes  in  which  radiative  interaction  from  base 
surfaces  could  be  ignored.  Herring  (10)  extended  the  same 
problem  to  include  specular  reflecting  surfaces. 

In  a  later  paper  (11)  Sparrow  and  Eckert  applied  the 
generalized  equations  to  a  configuration  consisting  of  a 
longitudinal  fin  connecting  two,  isothermal  cylindrical 
surfaces.  They  treated  the  surfaces  as  black  and  included 
all  radiative  interaction  effects.  An  iterative  Runga-Kutta 
technique  was  used  to  solve  the  governing  integral  differ¬ 
ential  equations.  Sarabia  (12),  by  application  of  an 
iterative  finite  difference  technique,  extended  the  work  of 
Sparrow  and  Eckert  to  include  grey  surfaces. 

Recently  Frost  and  Eraslan  (13)  published  a  study  of 
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radiation  from  arrays  of  longitudinal,  rectangular  profile 
fins  in  which  the  effects  of  mutual  irradiation  were 
considered.  They  solved  the  defining  non-linear  integro- 
differential  equation  through  an  iterative  technique  using 
the  variational  method  of  Galerkin.  Their  study  indicated 
that  interaction  effects  did  indeed  have  an  important 
influence  on  the  equilibrium  temperature  distribution  along 
the  fin.  Frost  and  Eraslan  also  concluded  that  the  Galerkin 
method  of  solving  the  differential  equations  had  a  distinct 
computational  advantage  over  the  heretofore  employed 
numerical  methods . 

The  study  of  the  combined  effects  of  convection  and 
radiation  from  single  fins  has  also  been  the  subject  of 
recent  investigations  (14,  15).  Frost  and  Eraslan  (16)  have 
extended  the  theory  to  arrays  of  longitudinal  rectangular 
fins  in  which  the  effects  of  radiative  interactions  under 
the  combined  modes  of  radiation  and  convection  heat  transfer 
are  considered.  Their  study  indicated  that  the  efficiency 
of  fins  under  the  combined  modes  of  heat  transfer  was  less 
than  for  the  pure  convection  case. 

The  present  analysis  extends  the  work  of  Frost  and 
Eraslan  to  cover  longitudinal  fins  of  arbitrary  cross- 
section  subject  to  the  combined  modes  of  heat  transfer.  The 
method  of  solution  is  formulated  to  include  the  effects  of 
non-uniform  surface  film  coefficients  and  surface  emissiv- 
ities  and  the  effects  of  mutual  irradiation.  However, 
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because  of  the  vast  amount  of  data  generated  by  an  exhaus¬ 
tive  study  of  this  type  of  fin  heat  transfer  phenomena,  the 
results  generated  by  the  method  of  solution  are  reported 
only  for  a  selected  group  of  fin  configurations  subject 
solely  to  radiative  heat  transfer. 
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CHAPTER  II 

FORMULATION  OF  THE  GENERAL  PROBLEM 

I.  GEOMETRICAL  CONSIDERATIONS  AND  ASSUMPTIONS 

A  typical  array  of  longitudinal  fins  of  arbitrary 
cross-section  projecting  from  a  plane  base  surface  is  shown 
in  Figure  1.  The  x-axis  is  chosen  to  run  parallel  to  the 
axis  of  the  fin  and  perpendicular  to  the  base  surface. 

The  analysis  is  predicated  on  the  following  assump¬ 
tions  : 

1.  Heat  conduction  along  the  fin  is  one-dimensional 
and  occurs  in  the  x-direction. 

2.  Each  fin,  while  of  arbitrary  cross-section, 
possesses  a  geometry  that  is  similar  to  all  other  fins  in 
the  array. 

3.  The  fins  are  immersed  in  a  non-participating 
radiation  medium  and  have  grey  surfaces. 

4.  Fin  surface  geometry  is  defined  as  a  function  of 
the  x-coordinate . 

5.  The  gas  medium  temperatures  are  known  functions 
of  the  x-coordinate. 

6.  Surface  film  coefficients,  thermal  conductivity, 
and  surface  emissivity  of  the  fin  are  prescribed  functions 
of  the  x-coordinate. 
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Fin  Cluster 


Figure  1.  Geometrical  conception  of  arbitrary  fin  cluster. 
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7.  The  heat  transfer  process  occurs  at  a  steady 

rate . 

8.  The  base  surface  exists  at  a  constant  tempera¬ 
ture  ,  T  . 

o 

9.  The  problem  is  formulated  for  a  fin  of  unit 

depth. 


II.  DERIVATION  OF  THE  BASIC 
DIFFERENTIAL  EQUATION 


The  differential  equation  of  the  phenomenon  is  formed 
in  the  classical  manner  by  making  an  energy  balance  on  a 
typical  differential  volume  element  as  shown  in  Figure  1. 

The  various  modes  of  energy  transport  considered  are 
internal  conduction  and  radiative  and  convective  transfer 
from  the  fin  surface.  Performing  the  energy  balance  pro¬ 
duces  the  following  form  of  the  general  equation: 


-K  (x)  A  (x) 


{■ 


+  <  K  (x)  A  (x) 


dT  (x) 
dx 


d^K(x)A(x)^y)-  dxj 


dx 


=  hx  (x)  A^  (x)  [T  (x)  -  Tgl(x)]  +  h2  (x)A2  (x)  [T  (x)  -  Tg2  (x)  ] 


+  qRl(x)A1(x)  +  qR2(x)A2(x) 


(1) 


where  K(x),  h1 (x) ,  and  h2  (x)  represent  the  fin  thermal  con¬ 
ductivity  and  surface  film  coefficients  of  the  corresponding 
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surfaces  respectively.  Expressing  the  area  terms  as  func¬ 
tions  of  the  fin  geometry,  simplifying,  and  collecting 
terms,  leads  to  Equation  2.  It  should  be  noted  that  fin 

J  A 

surface  slopes,  1  +  (dY/dx)  ,  are  included  in  Equation  2. 


Associated  with  Equation  2  are  boundary  conditions 
as  follows: 

1.  At  x  =  0 ;  T(0)  =  T  . 

o 

2.  At  x  =  L;  -K(L)^i^  -  qR*  +  h*[T(L)  -  T  (L)  ]  .  (3) 

y 

III.  CALCULATION  OP  THE  SURFACE  RADIATION  TERMS 
A  closer  look  will  now  be  taken  at  the  surface 
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radiation  terms  qRl(x)  and  qR2  (x)  .which  appear  in  Equations 
1  and  2.  The  radiative  energy  flux  emanating  from  the  fin 
surfaces  can  be  computed  by  treating  the  volume  between  the 
fins  as  a  cavity  and  applying  a  method  described  by  Sparrow 
and  Cess  (17)  for  analyzing  radiation  within  enclosures. 

The  method  of  Sparrow  and  Cess  assumes  that  the  walls 
of  the  enclosure  can  be  broken  up  into  (k)  isothermal  grey, 
opaque,  diffuse  surfaces.  The  radiosity  (B) ,  defined  as  the 
rate  of  energy  per  unit  surface  area  streaming  away  from  a 
given  surface,  is  composed  of  two  components;  the  emitted 

4 

energy  (eaT  ) ,  and  the  reflected  portion  of  the  incident 
energy.  The  energy  incident  on  unit  area  of  surface  (i)  is 
denoted  as  H^.  Figure  2  illustrates  the  radiosity  concept. 

Working  with  total  radiation  properties,  it  can  be 
shown  that  the  reflectivity  of  any  grey  surface  can  be 
represented  as  p  =  (1  -  e),  where  p  is  the  reflectivity  and 
e  is  the  emissivity.  Hence,  for  a  surface  denoted  with  the 
subscript  (i)  the  expression  for  radiosity  becomes  as  indi¬ 
cated  by  Equation  4. 

B.  -  e.aT. 4  +  (1  -  e. )  H.  (4) 

11111  '  ' 

Figure  3  illustrates  a  typical  surface  configuration 
in  the  enclosure  formed  by  adjacent  longitudinal  fins.  The 
i'th  surface  element  is  shown  receiving  radiation  from  other 
surface  elements,  from  the  base  region,  and  from  the  gas 


9 


AEDC-TR-69-88 


Radiosity 


Surface  i 


Figure  2.  The  radiosity  concept. 
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outside  the  cavity.  The  term  IT  in  Equation  4  is  seen  to  be 
the  sum  of  all  radiant  energy  arriving  at  the  i'th  surface 
element.  Denoting  any  other  surface  element  as  (j)  and 
noting  that  the  radiant  energy  which  leaves  (j)  is  in  fact 
the  radiosity  of  surface  (j)  leads  to  the  conclusion  that 
the  radiosity  of  any  surface  element  is  dependent  upon  the 
radiosity  of  all  other  surface  elements. 

The  radiation  energy  which  leaves  surface  element  (j) 
and  which  impinges  on  element  (i)  is  given  by  Equation  5. 


(Q„)  =  B.F . . 

R  j-i  ^  H 


(5) 


F is  the  shape  factor  which  relates  the  per  cent  of 
energy  leaving  surface  (j)  that  strikes  surface  (i) .  Summing 
the  radiative  contribution  of  all  surface  elements  to  the 
irradiation  of  surface  (i)  leads  to  Equation- 6. 


i  * 

h  =  ±  I  b.p.. 
1  Ai  j=i  3  3" 


.A.  =  T  B.F.  . 

11  jii  a  ^ 


(6) 


Substitution  of  Equation  6  into  Equation  4  and  solving 
for  the  temperature  term  produces  Equation  7 . 


.  B. 

CT . 4  =  — 
i  e . 

i 


(1  -  e. )  k  k 

- —  y  b.f.  ,  =  y 

1  1=1  J  J  1=1 


“  (1  -  ei)pi 


.  B. 
1  1 


£i 


(7)  J 


The  term  is  the  Kronecker  delta  which  takes  on  the  value 
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of  one  when  j  =  i  and  which  becomes  zero  when  j  ^  i. 

Equation  7  now  describes  a  system  of  (k)  linear 
algebraic  equations  relating  the  element  temperatures  and 
element  radiosities.  Expressed  in  matrix  notation,  the 
system  of  equations  is  represented  by  Equation  8. 

[aTi4]  =  [n±J]  [Bj]  (8) 

The  coefficient  matrix  elements  JL  ^  are  given  by 

=  [6..  -  (1  -  E.JP..1  /  e.  .  Inversion  of  the  coefficient 
matrix  gives  1  =  Ix^j).  By  employing  the  inverse 

matrix  a  new  set  of  linear  algebraic  equations  can  be  formed 
as  represented  by  Equation  9. 

[Bi]  =  txij]  toTj4]  (9) 

Equation  9  now  determines  the  radiosity  of  each  sur¬ 
face  element  as  a  function  of  the  temperatures  of  all  the 
surface  elements .  It  should  be  noted  that  the  open  end  of 
the  fin  enclosure  is  considered  a  hypothetical  surface  of 
unit  absorptivity  which  exists  at  the  effective  temperature 
of  the  gas . 

The  radiation  energy  flux  leaving  surface  element  (i) 
is  given  by  Equation  10. 

(qR}  =  e^T.4  "  (10) 
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Substitution  into  Equation  10  for  the  value  of  the 
incident  radiation  as  implied  by  Equation  4  gives  Equation 
11. 


(qR) 


i 


(11) 


Replacing  the  B^  term  in  Equation  11  with  Equation  9 
and  recalling  that  for  grey  surfaces  ot  =  e  gives  Equation  12, 
the  final  form  of  the  heat  flux  relationship. 


XijIaTj 


k 

-  I 

j=l 


Ai.cT. 


(12) 


The  term  A^..  in  Equation  12  is  seen  to  be  a  function  only  of 
the  particular  geometry  of  the  enclosure  and  the  surface 
properties  of  the  walls.  The  significant  aspect  of  Equation 
12  is  that  the  heat  flux  from  a  surface  element  can  be 
expressed  in  terms  of  the  discrete  distribution  of  tempera¬ 
ture  among  the  elements  which  make  up  the  fin  and  base 
surfaces . 


-  IV.  NON-DIMENSIONALIZATION  OF  THE 
BASIC  DIFFERENTIAL  EQUATION 

It  now  becomes  convenient  to  non-dimensionalize  the 
basic  differential  equation  (Equation  2)  and  its  associated 
boundary  conditions  (Equation  3).  Letting  the  parameters 
Kc,  and  hc  represent  characteristic  values  of  thermal 
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conductivity  and  heat  transfer  coefficient  respectively,  and 
taking  L  as  the  fin  length,  the  system  is  non-dimensionalized 
with  the  following  dimensionless  variables: 


v  =  £ 
X  L 


Y,  (x) 

vx>  -  -V- 


V2tx) 


Y2  (X)  =  L 


ic(x)  = 

c 


hl<X)  =  -jj 


h1(x) 


h9(x) 

Vx>  ■  -T— 


e  (x)  = 


T (x)  -  T, 

T  -  T 

geo  o 


T  ,  (X)  -  T 

9gl(X)  '  -  T 


goo  o 


T  _(X)  -  T 

eg2(X)  =  j  -  T 
*  goo  o 


h*  (i)  =  **£L 


qR1(x)  L 

qri(x)  =  iOf  :r~r  ) 

c  o  g°° 


qR9(x)  L 

°R2  (X)  =  K  (T  -  T  ) 
c  o  g°° 


S  = 


(13) 


Substitution  of  Equation  13  into  Equations  2  and  3 
gives  the  following  form  of  the  basic  differential  equation: 


f2(X) 


+  f,(x) 
dX^  A 


af  ’  fo(x)e<x> 


-Igc(x)  +  gR  (X)  ] 

(14) 


The  associated  boundary  conditions  upon  non- 
dimensionalization  become  as  shown  in  Equations  15. 
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Boundary  conditions: 

1.  0  (0)  =  0. 

2  d8(l)  -  Qp*  ^ 

-  ^kDT  ■  (15) 

The  functions  f 2 (X) ,  (X) ,  fQ(X),  gc (X) ,  and  gR(X) 

are  defined  by  the  following  equations: 


f2(X)  =  K  (X)  [Y1(X)  -  Y2(X)  ] 


(16) 


The  parameters  Nc  and  6*  are  defined  as  follows: 
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N 


C 


h  h*L 

A*  -  _£ _ _ 

P  KCK(1) 


The  details  of  transforming  Equations  2  and  3  into 
dimensionless  form  are  given  in  Appendix  A. 

Two  dimensionless  parameters,  the  radiation  number 
<Nr)  ,  and  the  temperature  ratio  Tgoo/To  are  now  introduced  by 
consideration  of  Equation  12.  Putting  Equation  12  into 
dimensionless  form  gives: 


(Q  )  =  N  l  A.  , 

R  i  j=l  3 


fT  1 

_i£  -  1 

To 

0  (X)  +1 

m 

w 

> 

(17) 


The  radiation  number  is  defined  by  N_  =  aT  L  /  K  (T  -  T  ) 

R  O  'co  goo' 


V.  TRANSFORMATION  OF  THE  NON-DIMENSIONAL  SYSTEM 


Inspection  of  the  second  dimensionless  boundary 
condition  (Equation  15)  indicates  that  this  boundary  condi¬ 
tion  has  a  non-zero  right-hand  side  and  hence  is  inhomoge¬ 
neous.  Application  of  the  B.  G.  Galerkin  technique  requires 
that  this  boundary  condition  be  homogeneous.  Hence,  it 
becomes  necessary  to  transform  the  dimensionless  system. 

Let  the  transformation  variable  [iMX)  ]  be  defined  by 
Equation  18. 


(X)  =  0(X) 


6*eg(l) 


K  (1) 


X) 


(18) 


Substitution  of  Equation  18  into  Equations  14  and  15 
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yields  the  following  transformed  system: 

f2(X)~"d'r'  +  £1  (x)^Sk"'  “  fo(X)*(X) 

-  |^*eg(i) 

-  [gc(X)  +  gR (X) ]  (19) 

The  associated  transformed  boundary  conditions  are  given  by 
Equation  20. 

1.  iMO)  =  0. 

2.  +  8*4>  (1)  =  0.  (20) 

Equations  19  and  20  now  define  the  mathematical 
system  which  will  be  solved  through  an  iterative  application 
of  the  B.  G.  Galerkin  method. 


3b 


X)fo(X)  - 


(2X-l)f1(X)  - 


2f 2 (X) 
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CHAPTER  III 

SOLUTION  OF  THE  DEFINING  EQUATION  BY 
THE  B.  G.  GALERKIN  METHOD 

I.  INTRODUCTION  OF  THE  GALERKIN  METHOD 

The  solution  of  the  mathematical  system  is  initiated 
by  assuming  that  Equation  19  has  a  solution  of  the  form 
indicated  by  Equation  21. 

n 

ijj  (X)  =  a..*.  (X)  +  l  a. 4>  <  (X)  (21) 

A  i=2  1  1 

where  <j>1,  and  <J> i  are  assumed  functions  of  X. 

Equation  21  must  satisfy  the  boundary  conditions  of 
the  problem  (Equation  19) .  Hence,  Equations  22  must  be 
satisfied. 

n 

<M0)  =  a. 4>  (0)  +  l  a.  4> .  (0)  =  0 
A  .  i=2  1  1 

- 3  +  3*^(1)  =  axta^d)  +  I  a.<j>.(l)] 

la  “  2 

n 

+  B*[a.<J>1  (X)  +  [  a. (1)] 

A  A  i=2  1  1 

=  0  (22) 
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A  form  for  the  functions  ^(X)  and  (X)  is  chosen  as 
indicated  by  Equation  23: 

<J»1(X)  =  (X  -  Y*)  X2 

*1(X)  =  (1  -  X)2  X1"1  (23) 


where  y*  =  (3  +  8*)  /  (2  +  8*). 

Appendix  B  demonstrates  that  the  assumed  functions 
and  do  indeed  satisfy  the  boundary  conditions  of  the 
problem.  The  problem  remains  of  evaluating  the  unknown 
coefficients  a^  and  a^(i  =  2,n)  appearing  in. the  solution. 

Substitution  of  the  assumed  solution  into  Equation  19 
yields  Equation  24. 


a1^f2(X)(6X  -  2y*)  +  f1vX)(3X2  -  2y*X) 

>}+  j2ai{[1 


-  fo(X)  (X3  -  Y*X2) 


:  i  <i  +  1)X 


-  2i(i  -  1JX1"2  +  (i  -  1)  (i  -  2)Xi-3]  f  2  (X) 
+  [  (i  +  DX1  -  2iX1-1  +  (i  -  l)Xi"2]f1(X) 

-  IX1+1  -  2X1  +  X1-1]fo(X)J 


-  -  mr) 


[(X*  -  X)  f  q  (X)  -  (2X  -  Df^x) 


-  2f  2  (X)  ]  -  [qc(X)  +  qR  (X)  ] 


(24) 
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Equation  24  is  now  a  functional  relationship  in  terms 
of  the  unknown  constants  a^  and  a^.  Following  the  method 
detailed  in  Kantorovich  and  Krylov  (18),  the  a's  can  be 
solved  for  as  the  dependent  variables  of  a  linearly  inde¬ 
pendent  system  of  algebraic  equations.  The  linear  system  is 
formed  by  multiplying  Equation  24  through,  in  turn,  by  each 
function  in  the  set  <Jk  (i  =  l,n)  and  integrating  the  result¬ 
ing  equations  between  the  limits  of  zero  to  one. 

II.  FORMULATION  OF  THE  ALGEBRAIC  SYSTEM 
OF  LINEAR  EQUATIONS 

Multiplying  Equation  24  through  by  <J>^  (X)  and  inte¬ 
grating  between  the  limits  of  zero  to  one  gives  Equation  25. 


a1[I2(l,l)  +  1^1,1)  -  I0(l,l)] 


n 


+  l  a.  [1,(1, i)  +  I.  (1 , i)  -  I  (1 , i>  3 
i=2  1  z  1  ° 


-{• 


V 


6*9  (1)  -  -  Ja(l)  -  J^l)]  -  Sll)  (25) 


The  terms  of  the  integrals  I,  J,  and  G  are  defined  as 
follows : 


I0(l,l>  =  I  fQ(X)(X3  -  y*X2)  2  dX 


■I 

^(1,1)  =  J  f± 


(X) (3X2  -  2y*X) (X3  -  y*X2)  dX 
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I2(l 


•i>  -r 


f2(X) (6X  -  2 y* ) (X3  -  Y*X2)  dx 


Io(1#i)  =  I  fo{X) (X1+1  "  2x1  +  xl  1) (x3  "  Y*X2)  dX 
''o 

=  I  f  1  (X)  [  (i  +  DX*"  -  2iXi-1  +  {i  -  l)Xi“2] 


*  (X3  -  Y*X2)  dx 


~r\ 


l2(l,i).  =  f  2  (X)  [i  (i  +  l)x1_i  -  2i(i  -  l)x3 


+  (i  -  1) (i, -  2)X1-3] (x3  -  Y*X2)  dx 


JQ(1)  =  (X2  -  X)  fo(X)(X3  -  Y*X2)  dx 


-r 

.{■ 
r2(l)  -J1 


(2X  -  1)  fx(X) (X3  -  Y*X2)  dX 


J2(l)  =  2f  2  (X)  (XJ  -  y*X^)  dX 


G{1) 


= J  gR(X) (X3  -  Y*X2)  dx  +  f  gc(X) (X3  -  Y*X2)  dX 

•'O  Jr\ 


Expressing  Equation  25  in  a  shorthand  notation  leads 


to  Equation  27: 


a,[Z(l,l)]  +  l  a.  [Z  (1  ,i)  J  =  ZH(1) 
i=2  1 
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The  elements  of  Equation  27  are  defined  as  follows: 

Z  (1,1)  =  [I2  (1,1)  +  I-l  (1,1)  -  I0(l,l)] 

Z(l,i)  =  [I2(l,i)  +  1^1,1)  -  IQ  (1 ,  i)  ] 

ZH(1)  =  |e*0g(l)  -  tao(l)  -  J2(l)  -  Jj.d)]  -  G  (1) 

(28) 

In  a  similar  manner  multiplying  Equation  24  through, 
in  turn,  by  each  of  the  functions  c|k  (i  -  2,n)  and  integrat¬ 
ing  between  the  limits  of  zero  to  one  produces  a  system  of 
algebraic  equations  represented  by  Equation  29. 


a1[I2(k,l)  +  I1(k,l)  -  l0(k,l)] 


n 


+  l  ai[I2(k,i)  +  I1(k,i)  -  IQ(k,i)] 

i  “  2 


■( 


6*0(1)  - 


r}[Jo 


^fiTpo(k)  -  J2(k>  -  Ji<k)I  -  G(k> 


(29) 


The  elements  of  Equation  29  are  defined  as  follows: 


V*'1) 


I^k,!) 


I2(k,l) 


■r 
-r- 

■r 


fQ(X)  (x3  -  Y*X2)  (xk+1  -  2Xk  +  xk_1)  dx 


(X) (3X2  -  2y*X) (Xk+1  -  2Xk  +  Xk_1)  dX 


f2(X) (6X  -  2y*) (Xk+1  -  2Xk  +  Xk-1)  dX 
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Vk 


-r 


fo(X)[X1+1  -  2X1  +  X1-1] 


(Xk+1  -  2Xk  +  Xk_1)  dX 


^(k,!)  =  j  fx 


I,  (k,i)  =  I  f ,  (X)  [  (i  +  ljX1  -  2iX1-1  +  U  -  DX1-2] 


(Xk+1  -  2Xk  +  Xk_1)  dX 


I2(k 


"  -r 


f  2  (X)  [  i  <  i  +  DX1-1  -  2i(i  -  1)X1-2 


+  (i  -  1)  (i  -  2) X1-3]  (Xk+1  -  2Xk  +  Xk_1)  dX 


Jo<k>  =|  (X2  -  X)f  (X)  (Xk+1  -  2Xk  +  Xk-1)  dX 


■r 

rl(k>  -  J  £x 


J, (k)  =  |  f, (X) (2X  -  1) (Xk+1  -  2Xk  +  Xk_1)  dX 


J2<k)  =  |  2f 2  (X)  (Xk+1  -  2Xk  +  Xk-1)  dX 


■r 

■r 

•r 


G (k)  =  |  gR(X)(Xk+1  -  2Xk  +  Xk_1)  dX 


+  |  gc(X) (Xk+1  -  2Xk  +  Xk_1)  dX 


(30) 


Expressing  Equation  29  in  a  shorthand  notation  leads 
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to  Equation  31: 

a.  [Z  (k,  1)  ]  +  l  a.  [Z(k,i)J  =  [ZH(k)]  (31) 

A  i=2  1 

The  elements  of  Equation  31  are  defined  as  follows: 


Z(k,l)  =  [I2(kfl)  +  Ix  (k,l)  -  IQ  (k, 1) ] 

Z  (k, i)  =  [  I2  (k ,  i )  +  I^kfi)  -  IQ(k,i)] 

zh  (k)  =  |e*0ga)  -  ^xr}[Jo(k)  "  J2{k)  "  Jl(k)]  “  G<k) 

(32) 

Equations  27  and  31  now  represent  a  system  of  n 
linear  algebraic  equations  containing  the  n  unknown  a^ 1 s . 
Expressed  in  matrix  notation  the  system  is  represented  by 
Equation  33 . 

n 

ZJl,m  am  =  ZHJ l  with  1  =  l,n  (33) 


which  upon  solving  for  the  a's  becomes 


n 


ZBm 

m=l 


(34) 
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III.  SOLUTION  OF  THE  ALGEBRAIC  SYSTEM 
OF  LINEAR  EQUATIONS 


The  Z&  matrix  elements  of  Equation  33  are  functions 

of  fin  geometry  and  the  assumed  form  of  the  (x)  functions. 

Hence,  once  the  fin  geometry  is  fixed,  the  Z„  elements  can 

*  ,m 

be  evaluated  and  the  inverse  computed. 

Equations  16,  26,  28,  and  32  indicate  that  the  [ZH] 

m 

elements  depend  upon  the  radiative  heat  transfer  from  the 
fin  which  is  in  turn  dictated  by  the  temperature  distribu¬ 
tion  along  the  fin.  Hence,  the  solution  of  the  linear 

system  of  equations  for  the  a's  requires  an  iterative 

process . 

A  simplified  outline  of  the  computation  steps 
necessary  to  achieve  a  solution  based  upon  a  selected  set 
of  fin  parameters  are  enumerated  as  follows: 

1.  Evaluate  the  matrix  elements  and  compute 

the  inverse  matrix  Z^m  for  the  selected  fin  geometry. 

2.  Assume  the  radiative  heat  transfer  from  the  fin 
is  zero  and  calculate  an  initial  temperature  distribution  in 
the  fin  based  on  convective  interactions  only. 

3.  Using  the  temperature  distribution  computed  in 
step  2,  evaluate  the  radiative  heat  transfer  from  Equations 
8  through  12  and  Equation  17. 

4.  Using  the  calculated  radiation  heat  transfer, 
evaluate  the  ZH  elements  and  calculate  a  new  temperature 
distribution  in  the  fin. 
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5.  Recalculate  radiative  heat  transfer  using  the 
temperature  distribution  from  step  4. 

6.  Repeat  steps  4  and  5  until  the  solution  converges. 
The  method  of  solution  presented  in  Chapters  II  and 

III  has  been  programmed  in  G  level  Fortran  on  the  IBM  360/50 
computer.  Appendix  C  is  a  Fortran  listing  of  the  computer 
program. 
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CHAPTER  IV 

APPLICATION  OF  THE  SOLUTION  TECHNIQUE 
TO  FINS  OF  SELECTED  SHAPE 

I.  DESCRIPTION  OF  SELECTED  FIN  SHAPES 
AND  CONDITIONS  OF  ANALYSIS 

The  foregoing  method  of  solution  was  applied  to  fins 
of  parabolic,  triangular,  and  inverse  parabolic  profile.  The 
dimensionless  area  of  each  profile  was  held  constant  at  0.2 
which  assures  the  same  volume  per  unit  depth  of  fin.  Thus  a 
comparison  of  the  reported  fin  performance  can  be  made  on  an 
equivalent  weight  per  fin  basis. 

The  selected  fins  were  analyzed  for  the  case  of  con¬ 
stant  thermal  conductivity,  negligible  convective  heat 
transfer,  and  a  ratio  of  gas  temperature  to  base  temperature, 
Tg«/To'  equal  to  zero*  These  latter  conditions  correspond 
closely  to  those  encountered  in  space.  Other  geometrical, 
environmental,  and  surface  property  parameters  varied,  as  are 
enumerated  below: 

1.  Non-dimensional  fin  spacing,  S,  was  assigned 
values  of  0.5,  1.0,  2.0,  4.0,  and  6.0. 

2.  The  radiation  number,  NR,  was  assigned  values  of 
0.1,  0.5,  1.0,  and  5.0. 

3.  Fin  surface  emissivity  and  base  emissivity  were 
set  equal  and  were  assigned  the  values  of  0.2,  0.5,  and  0.8. 
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The  total  number  of  cases  analyzed  was  180. 

II.  DESCRIPTION  OF  SELECTED 
MATHEMATICAL  PARAMETERS 

For  the  radiative  heat  transfer  calculations  the  fin 
surfaces  were  divided  into  forty  elements  and  the  base  area 
between  fins  into  five  elements.  The  number,  n,  of  approxi¬ 
mating  functions,  $,  used  in  the  Galerkin  solution  was  ten. 
Values  of  n  larger  than  ten  were  found  to  give  no  improvement 
to  the  first  four  significant  digits  in  the  solution  for  the 
temperature  distribution.  Figure  4  illustrates  the  effect 
of  variations  in  n  on  the  solution  for  temperature  at 
specific  locations  along  the  fin. 

To  facilitate  convergence  of  the  iterative  method,  0 
was  bounded  between  +2.0  and  -1.0.  The  physical  bounds  on  ‘ 

6,  +1.0  and  zero,  were  exceeded  in  the  solution  in  order  to 
stabilize  the  convergence;  since  as  the  solution  progressed, 
values  assigned  to  8  were  the  accumulated  average  of  all 
preceding  iterations. 

The  convergence  criterion  employed  required  tnat  the 
value  of  heat  flux  for  each  fin  element  should  change  by  no 
more  than  0.1  per  cent  during  successive  iterations.  The 
computer  was  programmed  to  run  for  thirty  iterations  or 
until  the  solution  converged.  Each  iteration  required 
approximately  three  seconds,  and  the  average  number  of 
iterations  required  per  solution  was  twenty. 
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Fin  Configuration 
Fin  Spacing,  S 
Surface  Emissivity,  e 
Radiation  Number,  NR 


Parabolic 


Figure  4.  Effects  of  n  on  the  solution  for  temperature 
distribution. 
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Cases  involving  high  values  of  surface  emissivity, 
high  Nr,  and  fin  spacing  greater  than  4.0  did  not  always 
converge  in  thirty  iterations.  These  cases  were  rerun  with 
fifty  iterations  and  in  all  cases  converged  to  an  acceptable 
solution.  Figure  5  shows  the  convergence  of  the  temperature 
distribution  for  a  typical  fin  configuration.  Bounding  of 
the  temperature  distribution  at  the  upper  bbund  is  indicated 
for  the  second  iteration  in  Figure  5. 
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Fin  Configuration 
Fin  Spacing,  S 
Surface  Emissivity,  c 
Radiation  Number,  Njj 


Inverse  Parabolic 
1.0 
0.2 
5.0 


Figure  5.  Convergence  of  successive  iterations. 
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CHAPTER  V 

PRESENTATION  AND  DISCUSSION  OF  RESULTS 

I .  TEMPERATURE  DISTRIBUTION 

The  non-dimensional  temperature  distribution,  T/Tq, 
for  each  combination  of  S,  NR,  and  e  are  shown  in  Figures  6, 
7,  and  8  for  the  parabolic,  triangular,  and  inverse  para¬ 
bolic  fins,  respectively. 

The  effect  of  fin  spacing  on  temperature  distribution 
is  seen  to  depend  upon  surface  emissivity.  For  low  surface 
emissivity,  e  =  0.2,  heating  of  the  fin  by  radiation  from 
the  base  is  appreciably  reduced  because  of  the  high  reflec¬ 
tivity  of  the  surface.  As  surface  emissivity  increases, 
fin  spacing  is  seen  to  have  a  greater  influence  on  temper¬ 
ature  distribution.  Physically,  this  effect  is  attributable 
to  increased  absorption  of  base  radiation  by  the  fin.  As 
spacing  increases,  the  fin  can  see  larger  areas  of  the  base, 
and  hence,  absorbs  more  of  the  radiation  from  the  base.  The 
overall  temperature  drop  along  the  fin  becomes  less  as 
spacing  increases;  since  as  radiation  energy  input  increases, 
the  conduction  of  energy  along  the  fin  must  decrease. 

Because  of  the  much  larger  shape  factor  between  the 
base  and  the  inverse  parabolic  fin,  a  greater  influence  of 
fin  spacing  is  observed.  It  is  interesting  to  note  that  for 
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Figure  6.  Dimensionless  temperature  distribution  in  a 
parabolic  fin. 
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Figure  6  (continued) 


Figure  6  (continued) 
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Figure  7.  Dimensionless  temperature  distribution  in  a 
triangular  fin. 
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Figure  8  (continued) 
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high  values  of  S,  NR,  and  e  the  temperature  gradient 
actually  reverses  near  the  two-thirds  point  on  the  fin. 

This  indicates  that  the  tip  region  is  heated  by  radiation 
from  the  base  at  a  rate  faster  than  it  is  cooled  by  reemis¬ 
sion  to  the  gas,  and  thus  heat  is  conducted  in  the  negative 
direction  along  the  fin  to  a  region  where  the  cooling  rate 
is  higher. 

Two  types  of  heat  transfer  phenomena,  conduction  and 
radiation,  occur  in  the  radiation  fin.  The  radiation  number, 
Nr,  being  a  ratio  of  radiation  to  conduction  effects,  is  a 
measure  of  the  dominating  mode  of  heat  transfer.  For  Nn 
approaching  zero,  the  heat  transfer  process  becomes  radia¬ 
tion  limited  since  the  fin  can  conduct  more  heat  internally 
along  the  fin  than  can  be  radiated  away.  Conversely  for 
Nr  greater  than  one,  the  process  becomes  conduction  limited. 
These  limits  are  clearly  demonstrated  by  Figures  6,  7,  and  8 
where  temperature  distributions  for  low  N  are  not  appre- 
ciably  affected  by  either  surface  emissivity  or  fin  spacing; 
while  for  large  NR,  a  strong  influence  on  temperature 
distribution  is  exerted  by  both  surface  emissivity  and  fin 
spacing. 


II.  FIN  EFFICIENCY 

So  that  the  performance  of  arrays  of  different  fin 
configurations  might  be  studied  and  compared,  it  becomes 
necessary  to  define  certain  fin  performance  parameters. 
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One  such  parameter  is  the  fin  efficiency,  n,  defined 
as  the  ratio  of  the  combined  convective  and  radiative  heat 
transfer  from  a  single  fin  to  the  hypothetical  combined 
convective  and  radiative  heat  transfer  from  a  fin  of  similar 
geometry  with  black  surfaces  and  with  infinite  thermal 
conductivity.  Expressed  in  terms  of  the  dimensionless 
variables,  the  fin  efficiency  becomes 


n 


f  {qri(x)  +  Ncu  "  Vx)1}  fT 


i  + 

’dYj^X) 

dX  i 

•I 


'R2 


if 

dY_ (X) 

+  Ncti  -  e2(x)]}  |i  + 

ax  J 

dX 


dX 


(35) 


Fin  efficiencies  are  calculated  for  each  fin  config¬ 
uration  and  presented  in  Figures  9,  10,  and  11  for  selected 
values  of  S,  NR,  and  e.  The  fin  efficiency  is  seen  to  vary 
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directly  with  surface  emissivity,  inversely  with  N_,  and  to 
increase  with  S,  asymptotically  approaching  some  maximum  and 
limiting  value.  The  limiting  value  corresponds  to  the 
efficiency  of  a  single  fin  situated  on  a  base  of  infinite 
dimensions,  and  free  of  interactions  with  adjacent  fins. 

Limiting  values  of  fin  efficiency. occur  at  a  spacing 
of  about  2.0  for  a  NR  of  5.0,  and  at  a  spacing  of  about  4.0 
for  a  Nr  of  0.1.  An  explanation  for  the  maximum  and  limit¬ 
ing  values  of  fin  efficiency  is  found  by  considering  that  as 
spacing  increases,  the  fins  interact  less  with  adjacent  fins 
and  more  with  the  warmer  base  region.  Increased  inter¬ 
actions  with  the  base  cause  the  overall  fin  temperature  to 
rise  because  of  higher  rates  of  radiation  energy  input  from 
the  base.  Further  spacing  cannot  influence  the  fin  tempera¬ 
ture  once  the  point  is  reached  where  adjacent  fins  no  longer 
interact,  since  the  fin  is  now  fully  interacting  with  the 
base.  At  this  point,  the  fin  has  a  maximum  temperature  and 
is  emitting  energy  at  the  highest  possible  rate  and,  hence, 
has  the  highest  possible  efficiency. 

Fin  efficiency  can  lead  to  serious  misinterpretation 
of  fin  performance  since  it  assumes  that  all  energy  leaving 
the  fin  also  leaves  the  enclosure  formed  by  the  adjacent 
fins  and  the  base.  When  interactions  are  present,  this 
assumption  is  false  and  can  lead  to  contradicting  conclu¬ 
sions  regarding  optimum  fin  configurations. 

A  more  valid  parameter  for  measuring  the  performance 
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of  arrays  of  longitudinal  fins  is  the  apparent  emittance. 

III.  APPARENT  EMITTANCE 

A  measure  of  the  performance  of  a  fin  cluster  is  the 
apparent  emittance,  eA<  Apparent  emittance  is  defined  as 
the  ratio  of  the  actual  rate  of  radiation  energy  efflux  from 
the  mouth  of  the  cavity  between  two  adjacent  fins,  to  the 
rate  of  radiation  energy  efflux  from  a  black  surface  at  the 
uniform  temperature  condition  of  the  wall,  Tq,  and  having  an 
area  equal  to  that  of  the  mouth  of  the  cavity.  Expressed  in 
terms  of  the  dimensionless  variables,  the  apparent  emittance 
becomes : 


Nr[1  -  (Tgoo/TQ)4]  {S  -  [Y1(0)  -  Y2(0)]}j 

(36) 


Apparent  emittances  were  calculated  for  each  fin 
configuration  and  presented  in  Figures  12,  13,  and  14  for 
selected  values  of  S,  N_,  and  e. 

X\ 

As  would  be  expected,  the  apparent  emittance,  eA,  is 
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strongly  influenced  by  surface  emissivity.  Higher  values  of 
surface  emissivity  are  seen  to  generally  yield  the  highest 
value  of  apparent  emittance  for  any  particular  value  of  S 
and  NR.  However,  a  reversal  of  the  cavity  effect,  i.e., 
ea  >  e,  for  high  emissivities ,  also  noted  by  Reference  (13) 
for  rectangular  fins,  exists  for  the  fin  configurations 
investigated  in  this  study.  For  e  =  0.2,  eA  is  generally 
greater  than  0.2.  However,  as  c  increases  to  0.5  and  0.8, 
eA,  in  most  cases,  becomes  less  than  e.  The  explanation  for 
this  trend  lies  in  the  fact  that  for  low  values  of  e,  the 
base  radiation  is  reflected  from  the  fins  and  augments  the 
radiation  to  the  external  environment,  while  for  high  values 
of  e,  the  base  radiation  is  absorbed  by  the  fin  surfaces  and 
is  reemitted  at  a  lower  temperature  resulting  from  the 
conduction  loss  along  the  fin,  thus,  limiting  the  base 
contribution  to  the  overall  heat  transfer  rate. 

The  influence  of  NR  upon  eA  for  a  fixed  value  of  fin 
spacing  is  also  apparent  in  Figures  12,  13,  and  14.  The 
apparent  emittance  is  seen  to  decrease  as  NR  increases 
because  either  the  amount  of  base  radiation  absorbed  by  the 
fins  becomes  increasingly  larger  or  internal  conduction 
along  the  fin  is  decreased. 

A  comparison  of  ea  for  the  various  fin  configurations 
is  presented  in  Figure  15.  Curves  of  ea  are  presented  for 
Nr  values  of  0.1  and  5.0,  and  for  e  values  of  0.2,  0.5,  and 
0.8.  Data  from  Reference  (13),  for  a  rectangular  fin  of 
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Figure  15.  Apparent  emittance  of  clusters  of  fins  of 
selected  shapes. 
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equal  non-dimensional  volume  per  unit  length  and  fin  depth, 
are  presented  for  comparison  purposes.  Surprisingly  little 
difference  is  observed  between  the  parabolic,  triangular, 
and  inverse  parabolic  profiles.  The  rectangular  fin  is  seen 
to  have  the  highest  value  of  eA  for  any  selected  set  of  S, 
Nr,  and  e.  The  data  from  Reference  (13)  for  the  rectangular 
fin  have  not  been  verified  using  the  computer  program 
developed  for  this  study. 
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CHAPTER  VI 

CONCLUSIONS 

The  application,  by  Reference  (13),  of  an  iterative 
technique  to  the  B.  G.  Galerkin  method  of  treating  the 
defining  integro-dif ferential  equation  constitutes  a  highly 
successful  approach  to  the  solution  of  this  class  of 
problems.  The  current  extension  of  the  method  to  include 
fins  of  arbitrary  profile  was  equally  successful. 

Application  of  the  method  of  solution  shows  that 
radiative  interactions  can  strongly  influence  the  temperature 
distribution  and  heat  transfer  characteristics  of  clusters 
of  longitudinal  fins. 

Analysis  of  results  indicates  that  only  under  certain 

circumstances  can  fins  enhance  the  radiation  heat  transfer 

from  the  base  region.  For  e  k  0.8  no  advantage  can  be 

gained  by  using  fins.  Similarly,  for  e  =  0.5  and  N  >  0.5 

R 

no  advantages  accrue  from  using  fins,  however,  for  N  =0.1 

R 

an  improvement  in  heat  transfer  characteristics  occur  when 
fins  are  used.  In  the  case  of  surfaces  with  e  £  0.2,  fins 
will  improve  the  heat  transfer  from  the  base  region  for  all 
values  of  NR  £  1.0. 

When  conditions  indicate  that  fins  will  improve  the 
heat  transfer  from  the  base  region,  investigation  reveals 
that  for  maximum  benefit  the  fin  spacing  should  be  as  close 
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as  possible.  This  conclusion  Contradicts  the  conclusion  one 
might  draw  from  inspection  of  the  fin  efficiency  curves, 
thus,  illustrating  why  fin  efficiency  is  not  a  valid 
parameter  for  fin  comparisons  when  radiative  interactions 
are  present. 

The  rectangular  fin  profile  examined  by  Reference  (13) 
appears  to  be  superior  to  either  the  parabolic,  triangular 
of  the  inverse  parabolic  profile  when  compared  on  the  basis 
of  equal  volume  per  unit  non-dimensional  depth  of  fin. 
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APPENDIX  A 

NON-DIMENSIONALIZATION  OF  THE  BASIC 
DIFFERENTIAL  EQUATION 

Frequently  it  becomes  advantageous  to  express  the 
solution  of  a  differential  equation  in  dimensionless  form. 
There  are  at  least  two  reasons  for  non-dimensionalizing  an 
analysis  as  enumerated: 

1.  A  study  conducted  under  a  particular  set  of 
parametric  conditions  can  be  generalized  to  apply 
to  many  sets  of  parameters  by  non-dimensionalizing 
the  solution. 

2.  The  effects  on  the  solution  of  particular 
parametric  groups  can  be  isolated  and  more  easily 
studied  if  the  system  is  non-dimensionalized. 

Equations  2  and  3,  representing  the  dimensional 
differential  equation  and  its  associated  boundary  conditions 
respectively/  are  rewritten  from  the  text  and  presented  as 
Equations  A-l  and  A- 2. 
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(A- 2) 


The  dimensional  variables  of  Equation  A-l  are  non- 
dimensionalized  by  introducing  the  following  set  of 
dimensionless  variables: 
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Performing  some  initial  manipulation  of  variables  in 
the  differential  terms  gives  the  following  results: 

T(x)  =  (Tgoo  -  Tq)  0  (X)  +  Tq 

dux)  Jny  -  v  e<*> + 10]  ax 

dx  <JX  dx 

but 


dx  _  1 

dx  ”  L 


hence 


dT(x)  _  (Tq«>  ~  Tp)  d6  (x) 
dx  L  dX 


d2T (x) 


(Tgco  ~  Tq>  d0  (X) 


dx 


(T 

ax  _  v  g 


dX 


dx 

dx 


-  V 


d2e  op 
dx2 


Substitution  of  the  non-dimensionalized  differential 
terms  together  with  the  terms  of  Equation  A-3  into  Equation 
A-l  yields: 


(T 


kCxJk^y^x)  -  y2(X)]l  —3. 
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r°°  O 

1 - 


d20 (X) 
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Factoring  out  the  constant  terms  that  were  introduced 
in  the  non-dimensionalization  process  gives  the  following 
equation: 
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Next  by  dividing  Equation  A- 5  through  by  the  factored 
constant  group  K  (T  -  T  )/L,  the  following  transformed 
equation  is  derived: 


+  [Y1(X)  -  Y2(X)] 


de  (x) 
dx 
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Equation  A-6  is  the  non-dimensional  form  of  the 
general  differential  equation  and  is  written  in  the  text  as 
Equation  15. 

The  first  boundary  condition  is  non-dimensionalized 
by  a  simple  substitution  of  the  dimensionless  variables  into 
Equation  A-2.  Equation  A-7  is  the  non-dimensional  expres¬ 
sion  for  the  first  boundary  condition. 

T  ( 0)  -  T  T  -  T 

At  X  =  0,  0  (0)  =  5 - -f-S-  =  ^  P  ■  =  0  (A-7) 

g«o  o  g°°  o 

Upon  substitution  of  the  dimensionless  variables  into 
the  second  boundary  condition,  the  following  relationship 
is  produced: 
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Now  by  letting  6*  =  hch* (l)L/K(l)Kc,  the  second 
boundary  condition  takes  on  its  final  dimensionless  form: 

+  6*e(D  =  6*0  (i)  -  %-&-  (a- 8 

dX  *  K  (1) 

Equations  A-7  and  A-8  are  represented  as  Equation  15 
in  the  text. 
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APPENDIX  B 

BOUNDARY  CONDITION  CHECK  OF  ASSUMED 
GALERKIN  FUNCTIONS 

It  is  demonstrated,  in  the  text,  how  the  assumed  form 
of  the  B.  G.  Galerkin  functions  are  applied  in  obtaining  a 
solution  of  the  differential  equation.  It  remains  to  show 
that  the  assumed  form  of  the  solution  does  satisfy  the 
boundary  conditions. 

The  assumed  solution  and  the  boundary  conditions  are 
rewritten  from  the  text  and  stated  here  as  Equations  B-l  and 
B-2  respectively. 

Assumed  Solution: 

n 

iMX)  =  a.^.U)  +  l  a.<t> .  (X)  (B-l) 

1  i=2  1  1 

Boundary  Conditions: 

1.  \^(0)  =  0  at  X  =  0 

2*  +  B*iH1)  =  0  at  X  =  1  (B-2) 

The  assumed  form  of  the  Galerkin  functions  is  repeated 
here  in  Equation  B-3. 
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4>1  (X)  =  (X  -  y*)  x2 
Y*  =  <3  +  B*)  /  (2  +  B*) 

^(X)  =  (1  -  X)2  Xi_1  (B-3) 

Substitution  of  the  assumed  solution  into  the  first 
boundary  condition  gives: 

n 

iM0)  =  a,  (0  -  Y*)  (0)2  +  l  a.  (1  -  0)2  (0)i_1  =  0 

i=2  1  (B-4) 

Next,  substitution  of  the  assumed  solution  into  the 
second  boundary  condition  gives: 

+  B*iMD  = 

dX 

+  B*  [a,  <f>1  (1)  +  l  a.  <)>.  (1)] 
i=2 

-  a^d)2  -  2  (1) y*] 
n 

+  .l  at[  (i  +  1)  (l)1  -  2i  (1) i_1  +  (i  -  l)  (l)1"2] 

i-2 

+  B*{a1[(l)3  -  (1)V]  +  I  a.[U>1+1  -  2d)1  +  (l)1"1] 

i=2 

■  axt3  -  2y*]  +  B*[a1(l  -  y*  )  ] 
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3^  3-2 


3  +  6*1 

[mrj 


+  3*  ax  1  - 


3  +  3* 
2  +  3* 


-  j-2  +\*)  [6  +  33*  -  6  -  23*  +  3*(2  +  3*  -  3  -  3*)]  =  0 

(B—  5 ) 

Thus  it  is  shown  how  the  assumed  form  of  the  Galerkin 
functions  satisfy  the  boundary  conditions. 
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APPENDIX  C 

COMPUTER  PROGRAM  LISTING 


C  GENERAL  SOLUTION  FOR  RADIATION  CONVECTION  FIN 

IMPLICIT  RFAi.*RIA-H,0-Z> 

DIMENSION  Z!30.30>,ZH!30>.ZH2I30).COLI)MN!30),SA!30)»CHII60,60), 
1CAPLAMI 60*60) *QUEOOT !  50)  »  QL AST  I  50 ) *C APPS I ( 30 )« 

20 IFF  150  >* EX ( 30) *THET  Al 30 ) • A! 30 ) * ZH1 I  301 » SHAPE (60*60) *TRAT! 30) * 
3TEMPAV ( 50 ) 

DIMENSION  F2XI51  ),F1X( 5l),F0X<51>  ,GXI 51) 

SORT IO)=DSQRT (0) 

ARS I P  >=DABS  t 0 ) 

Y1X|X  )=. 3* ( 1 • -X  >**2 
Y2X(X)=-.3*I1.-X)**2 
DYlDX!X)=-.6*< l.-XJ 
0Y20X  !  X  >  =  .6*11. -X) 

H1X(X )=L.+X-X 
H2XIX  )=1.+X-X 
SY1(0X1)=1.+(0X1*0X1) 

S Y2 ( 0X2 )  =  l. +10X2*0X2) 

F2X1! Y1,Y2)=AKX*!Y1-Y2)*AR 

F1X1! Y1,Y2,DX1,0X2)*!DKDX*!Y1-Y2>+AKX*IDX1-DX2) l*AR 
F0X1CH1,H2,S1.S2)=B*{ IH1*!S1**0.5))+!H2*IS2**0.5) ) )*AR 
GX1!H1*H2*S1*S2)=B*( ( H1*TG1X*! S 1**0 . 5 )) +1 H2*TG2X*! S2**0. 5 )) )*AR 
EPS2  <  X >=EFI N+X-X 
EPSHX).=  EFIN+X-X 
B=0.0 

WRITER. 1110) 

1110  FORMAT! ■  INVERSE  PARABOLIC  FIN  WITH  SLOPE  INCLHDEO. •) 

WRITE! 6*1120) 

1120  FORMAT!*  UNIFORM  FILM  COEFFICIENT  ,HIX=H2X=1.0  •) 

C  DEFINE  FIN  FUNCTIONS 

AR=0. 04/3 . 

AKX=1. 

AKO-AKX 

AKC=AKO 

DKDX=0. 

TG1X=1. 

TG2X=1. 

C  EVALUATION  OF  FIN  FUNCTIONS  AT  X=0. 

X=. 00000001 
0A0=Y1X ( X  >— Y2  X I  X ) 

F2XI1  )*>F2XHY1X!X),Y2X!X  )  >*0.5 
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FI  XU >=F1X1(Y1X(X),Y2X(X) , DY 1DX (X > .DY20X (X> 1*0.5 
F0X(1>=F0X1<H1X<X>,H2X(X>,SY1(DY1DXIX1  >,SY2(0Y2DX(X)  >  J*0.5 
GX( 1 (=GX1 IH1XIX  >,H2X( X  >  tSYXI DY1DXI X> ) »SY2(0Y20X( X ) ) >*0.5 
END  X=0. 

EVALUATION  OF  FIN  FUNCTIONS  AT  X=l. 

EVALUATION  OF  TIP  PARAMETERS 
X=. 99999999 
DAS=Y1X(X>-Y2X(X> 

NS=0.5*(H1XIX>+H2X(X>> 

BS* I HS/AKX ) *B 
TGS«*0.5*  CTG1X+TG2X) 

GMS=l3.+BS>/(2.+BS> 

F2XI51 >=F2X1 ( Y1X I X>.Y2X(X>)*0.5 

FIX  (51  >*F1  XI  mX(X»,Y2X(X>,DYlDX(X>,DY20X(X  >>*0.5 
FOX(51>=FOXI(H1X(X>,H2X(X>,SY1(DYIDX(X>  > , SY2 1 0Y20X ( X ) ) >*0 . 5 
GX ( 51 >=GX1(H1X(X>,H2X(X),SY1(0Y10X(X>  J, SY2(DY2DX( X) >  >*.5 
ENO  X«l. 

END  DEFINITION  OF  TIP  PARAMETERS 
SX=-1. 

00  1030  IX-2.50 

sx=-sx 

AX=l. +0.5*1 l.+SX> 

XI*1 X-l 
X=0.02*X  I 

F2X(  IX)  =  F2X1(Y1X(X),Y2X(X)»*AX 

=  IX( IX  »  =  F1X1 (YlX(X)fY2X(X)»DYlDX(X)tDY2DX(X))*AX 

FOXI IX)  =  F0X1(H1X(X),H2X(X> , SY1 ( DY 1DX < X >  > , SY2 ( DY20X ( X > ) >*AX 

GX( !X>=GX1(HIX(X) .H2X(X),SY1(0Y1DX(X>  ) , SY2< 0Y2DX { X > ) )*AX 

CONTINUE 

ENO  FUNCTION  STORAGE 

MAXIMUM  MATRIX  DIMENSIONS,  NC=20,  ND=21 

START  LOOP  FOR  EVALUATION  OF  MATRIX  PLFMENTS 

NHel 

NC=20 

NO*NR+NC 

S 1 2 1 1=0. 

SI111=0. 

SI011=0. 

sjoi=o. 

SJ11=0. 

SJ21=0. 

SSHF=0. 

SJS1-0. 
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C  START  MUMFRICAL  INTEGRATION  FOR  (1,1)  AND  (l,NQ) 

1)01031  IX=1  ,51 
X 1= IX-1 
X=0  »02*X I 

phx=x*x*x-gms*x*x 

TI211=(F2X( IX) )*I6.*X-2.*GMS)*P11X 
T!111=(F1X( IX) )*(3.*X*X-2.*GMS*X)*PllX 
TIOU  =  (FOX(  IX  >  )*(X*X*X-GMS*X*X  >*P11X 
T  J01= ( FOX  I  IX) )*(X*X-X)*P11X 
TJ11=(FIX(IX) >*(2.*X-1 ,)*P11X 
TJ21=(F?X( IX)  >*P]1X 
T  JS1  =  (GX ( IX) )*P11X 
TSHF=GX( IX) 

SI2ll=SI211+TI2ll 
Sllll-SIll 1+T 1 1 1 1 
S 101 1 =S 101 1+T 1 01 1 
SJ01=SJ01+TJOl 
S J 1 1  =  S Jll+T  J1 1 
SJ21=SJ21+TJ21 
SJS1=SJSI+TJS1 
SSHF=SSHF+TSHF 
1031  CONTINUE 

C  ENO  INTEGRATION  FOR  (1,1)  ANO  (l, NO) 

Z(1,1)=S 121 1+SI111-S 1011 

Z»( 1 )=RS*TGS*(SJ01-SJ11-2.*SJ21 )-SJSl 

ZH1(1)=  (SJ01  -  SJll-  2.*  SJ21  ) 

SHF=SSHF 
005000  1=2, NC 
ZI=1 

S121 1=0. 

Sill  1=0. 

SI01I=0. 

C  START  NUMERICAL  INTEGRATION  FOR  (1,1) 

001032  IX=l ,51 
XI=IX-1 
X=0.02*XI 
PI IX=X-GMS 

TI21I=IF2X(IX))*(( (ZI-1.)*(ZI-2.)*(X**< 1-1) ))-( (2.*ZI)*(ZI-l.)*(X 
1**1  ))  +  (  (ZI+1.-)*(7I  >*(X**(I+1) ) )  )*P1IX 

TI111=(FIX(IX))*(((ZI-1.)*(X**I ))-((2.*ZI )*(X**( I+l) ))♦( (ZI+1.)*( 
1X**(!+2)))I*P1IX 

TIOU=(FOX( IX>)*( (X**( I  +  l) )-(2.*(X«* ( 1+2) ) )+IX**{ I+3>))*PLIX 
S 121 I=SI21 1+T 121 1 
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SIllI=SUl  1+T1111 

sinu=sioii+Tioir 

in 32  CflftiT  I‘Vl)P 

C  FND  I  NT  t-RR  AT  1  (IN  POP  11,1) 

ZI1,I >=51211 +SI11I-S 1011 
SOOO  CONTINUE 

005010  <=?,nc 

7K=Ki 
S I 2K 1 =0. 

S I 1KI=0. 

Sin<i=o. 

sjok=o. 

SJ1<=0. 

SJ2K=0. 

SJSK=0. 

C  START  NUMERICAL  INTEGRATION  FOR  (K,l>  AND  IK, NO) 

001033  IX=1 ,51 
X 1= I X-l 
X=0.02*XI 

PK1X=<  X**IK-l )  )-I2.*(X*«‘K>>  +  (X**(K+l)  t 
TI2<1=(F2XI1X) >*(6.*X-?.*GMS)*PK1X 
TI1K1=IF1X( IX) )*I3.*X*X-2.*GMS*XJ*PKIX 
TI0K1  =  IF0XI1X>>*IX»'*3-GMS#X*X>*PKIX 
T  JOK  = I  FOX  I  I X ) ) *( X*X-X ) *PK IX 
TJ1K=(FIXIIX>>*(2.*X-1.)*PK1X 
TJ2K=|F2X( IX  >  >*PK1X 
TJS<=(GXI IX))*PK1X 
SI2K1=SI2K1+TI2K1 
SI1K1=SI1K1+TI1KI 
SIOK1=SIOK1+T!OK1 
SJOK=SJOK+TJOK 
SJ1K=SJ1K+TJ1K 
SJ2K=SJ2K+TJ2K 
SJSK=SJSK+TJSK 
1033  CONTINUE 

C  ENO  INTEGRATION  FOR  ■( K ,  1  >  AND  IK, NO) 
ZIK,1>=SI2K1+SI1K1-SI0K1 
ZH(K)=BS*TGS*ISJ0K-SJ1K-2,*SJ2K>-SJSK 
ZH1IK)  = I SJ0K-SJ1K-2 ,#5J2K  > 

5010  CONTINUE 

005020  K  =  2  ,NC 
ZK=K 

005030  1=2, NC 
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7.1  =  1 

S l 2< 1=0. 

SI1K 1=0. 

SI0KI-0. 

C  START  NUMERICAL  INTEGRATION  FOR  IK, I) 

DO 10 34  1X=2,51 

XI-IX-1 

X=0.02*XI 

C  IX  STARTS  AT  2  SINCE  FOR  IX=1  X=0  ANO  THUS  PKIX=0 
P*1X  =  <  <1.-X)*(1.-X)>MX**(  I+K-4H) 

TI2KI  =  (F2XIIX>'»*U(ZI-l.>*(ZI-2.)>-CC2.*Zn*IZI-l.)*X»  +  |(2I+l.I*ZI 

1*X*X  > )*PK I X 

T  I  IK  !  =  <F1X(  I X » >*UZ  1-1.  >-f 2.*Z1*XI+I<ZI+1.  )*X*X) )*X*PKIX 
T!OKt  =  (FOXUX>)*<l.-X)*(l.>X)*X*X*PKIX 
SI2KI=SI2KI+TI2KI 
SI1KI=SI1KI+TUKI 
SI0KI=SI0KI+T10K1 
1034  CONTINUE 

C  END  INTEGRATION  FOR  IK, I) 

Z(K,I>=SI2K1+SI1KI-SI0KI 
5030  CONTINUE 
5020  CONTINUE 

C  END  EVALUATION  OF  MATRIX  ELEMENTS 

C*********************************************** 

C*********«*******«***«*M******************* 

C  THE  FOLLOWING  ARE  CERTAIN  INOEXING  CONSTANTS  AND  SOME  INITIAL  VALUES 
1=20 

C  I  REPRESENTS  THE  NUMBER  OF  ELEMENTS  THE  FIN  IS  DIVIDED  INTO  FOR  RAD.INTEG. 
K  =  I 

C  N  REPRESENTS  THE  ORDER  OF  THE  GALERKIN  ASSUMPTION  ON  THE  SOLUTION 
N=10 

C  IB  REPRESENTS  THE  NUMBER  OF  ELEMENTS  INTO  WHICH  THE  BASE  IS  DIVIDED 
IB=5 
12  =  1+1 
13  =2*1 
14=13*1 
1 5= 1 3* IB 
14=15+1 

C  IEND  CORRESPONDS  TO  THE  TOTAL  ELEMENT  NO* AND  INCLUDES  BASE 
!FND=  15 
JSTOP= 1 6 
THET  A (K )=0.0 
OSTAR  =0.0 
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TRATC»0.0 

EPSGAS-.9999 


THE  FOLLOWING  IS  A  SEARCH  FOR  THE  SMALLEST  ELEMENT  OF  THE  Z  MATRIX 
AMINoi.o 
00  350  I Y"1 «N 
DO  350  KYb1,N 

AAbDABS(Z(IY,KY))-DABS(AMIN> 

IF  ( AA >  351,351*350 
351  AMINbOABSIZI IY,KY>> 

350  CONTINUE 


THE  FOLLOWING  IS  A  NORMALIZATION  OF  THE  Z  MATRIX  BY  DIVISION  BY  THE  MIN  ELEM. 
THE  MINIMUM  ELEMENT  IS  CALLEO  AMIN 
DO  180  N1=1,N 
DO  180  N2*l ,N 
ZCN1.N2)  =Z (N1 ,N2  >/AMIN 

180  CONTINUE 

C  ************************************* 

C  CALLING  THE  REVISEO  DOUBLE  PRECISION  INVERSION  ROUTINE"VERTDR" 

CALL  VERTDR  (N,30,Z,ICHECK> 

£********************«******•*****«********* 

00  181  N3»1,N 
DO  181  N4sl,N 
Z I N3»N4)*Z (N3*N4)/AMIN 

181  CONTINUE 

C*******************»****«*****************« 

IF  UCHECK  +1)  133,134,133 
134  WRITE  (6,610) 

610  FORMAT! •  ZIK,I>  MATRIX  IS  SINGULAR  -  SORRY  ABOUT  THAT  CHIEF* I 
133  CONTINUE 

C  ft************:************************ 

C*************************4************«***e 

C4****************************************** 

EFIN-.2 

C  INSERT  DO  HERE  FOR  PERMUTATION  ON  FIN  EMISSIVITY 
8ASE=1. 

C  INSERT  DO  HERE  FOR  PERMUTATION  ON  FIN  SPACING 

C  BASE  IS  THE  NONDIMENSIONAL  FIN  PITCH 

C  ALONG  IS  THE  RATIO  OF  NON-DIM. FIN  LENGTHIN  SUBROUTINE  FACTOR  TO  THF  NON- 
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C  DIMENSIONAL  FIN  LENGTH  IN  THE  MAI*  PROGRAM, 


ALONfi=l./ (HASP-. 6) 


C  RA5E1  IS  The  RASE  WIDTH  Of  THE  UPPER  PIN  SURFACE. 

BASE1=  .3*AL0MG 

C  BASE2  IS  THF  RASP  WIDTH  OF  THE  LOWER  fin  SURFACE. 

BASE2=  .3#ALONG 

(!♦  ♦***  +  #*#***#*  **  4*6*4  **  *$*#*$ 

C?****  **#**###*#**<1 

C******************’********** l**#*#*#****#**^*!),##**,,,,,,,^,  . 

C  THE  FOLLOWING  SECTION  TAKES  THE  SHAPE  FACTOR  MATRIX  'SHAPE' 
C  SUBROUTINE  FACTOR  AND  COMPUTES  THF  'CHI'  MATRIX 
CALL  FACTOR  (SHAPE  , I , I B .BASF  1 ,BASE2, ALONG } 

C  CALCULATION  OF  THE  'CHI'  MATRIX  ELFMENTS  FOLLOWS 
DO  171  JL  =  1 ,16 
00  160  I L  =  1 *  I 
CI  =  I 
C IL  =  IL 

EPS  11=  EPS1( (2.*CIL-1.  >*(  l./(2.*CI )  >  I 
IF  (IL-JU  161,162,161 

161  CHI  (  IL  ,  JL  >=  (-( 1  .-EPS  I  D/EPSll)*  SHAPEUL.JU 

GO  TO  160 

160  CONTIn'e1**  =  <1*“<  *  1,~EPSI1  >*SHAPE(  IL,  JL)  )  I/EPSI1 
00  163  IL=I2,I3 


164 

165 
163 


167 

168 
166 

169 

170 

171 


EPS  12=  EPS2(  (2.*(CIL-CI>-l.)*(l./(2.*cnn 

IF  C IL— JL )  164,165,164 

CHI(!L,JL>  =  (-(l.-EPSI2)/EPSI2>*SHAPE( IL, JL) 

GO  TO  163 

CHI (IL,JL  >  ■( l.-r( (l.-EPSI2>*SHAPE( Ib.JL)  IJ/FPSI2 

CONTINUE 

DO  166  IL=  14,15 
EPS3=EFIN 
EPSI3=  EPS3 
IF  (IL-JU  167,168,167 

CHI (IL, JL)  =  (-( l.-EPSI 3  >/EPSI 3 )*SHAPE ( IL , JL I 
GO  TO  166 

*  =  <1--H1.-EPSI3>*  SHAPE(  IL,JL)))/EPSI3 

CONTINUE 

IF  (I6-JL)  169,170,169 

CHI ( 16, JL )  =  (-( l.-EPSGAS I/EPS GAS)*  SHAPE(I6,JL) 

GO  TO  171 

rni!!fe’JL>  =  <!•-< I l.-EPSGAS»*SHAPE( 16, JL) II/EPSGAS 
CUN  1  I NUc 


FROM 
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C  calling  sif»wniiTINP  vphI'i 

CALL  VEPTO  < 16.60, CHI ,JCHECK ) 
i r  ur.Hcc<+i)  i7?.i n,  17? 

173  >JRI1F|6,630) 

630  FORMAT  I'  THE  CHI  MATRIX  IS  S I  NGULAR  •  THE  PROGRAM  WILL  CONTINUE') 
17?  CONTINUE 

Q ****************************************************** 
£*********************************************** 

C  EVALUATION  OF  THE  CAP  I. A  1*0 A  FUNCTIONS  USING  THE  CHI  MATRIX 
DO  100  M  =  1 ..!  STOP 

C  THE  SUBSCRIPT  M  IS  EQUIVALENT  TO  THF  SURSCRIPT  I  IN  THE  ANALYSIS 

CM=M 

IFIM-I)  103,103.101 

101  IF ( M-I 3)  104,104,102 

102  1FIK-I5)  105 ,103,106 

103  EPS  =  EPS1 ( (CM-.51/CI ) 

GO  TO  107 

104  EPS  =  EPS2I1ICM-CI  J-.5I/CI  ) 

GO  TO  107 

1.05  EPS  =  FPS3 
GO  TO  107 

106  EPS  =  FPSGAS 

107  EPSCflN  =  EPS/ll. -EPS) 

00  100  MM=1 , JSTOP 

C  THE  SUBSCRIPT  MM  IS  EQUIVALENT  TO  THE  SURSCRIPT  K  IN  THE  ANALYSIS 

CN=MM 

IF  1008, 109, 100R 

1008  CAPLAM|M,MM)=  FPSCON  *  ( -CHI ( M, MM) ) 

GO  TO  100 

109  C A PL AMI M,MM ) =  EPSCON  *  I l .-CHI  I M, MM ) > 

100  CONTINUE 

C*********************************************** 

C  INSERT  on  HERE  FOR  PERMUTATION  ON  RADIATION  NUMBER  IRN). 

RN=-5. 

C ************************* ** *** *********** ****** 

icniiNT  =  o 

£*«*** ****************************************** 

00  108  M=l, JSTOP 

0UE0nT|M>=0.0 

TEMPAV|M>=0.0 

108  CONTINUE 

C*********************************************** 

GO  TO  150 
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2000  CONTINUE 

no  1001  Msl , JSTOP 
SIPM=0.0 

r>f)  1101  MHO  1,1 

SllM  =  StlM+CAPL  AM(M,  MM)*M1.+  (ITRATC>  — 1.0)*THFTA(MM1  )**4) 

1101  CONTINUE 

on  110?  MM=  I ? , I  3 

SIIM=SUM+CAPLAM(M,MM)*( I l.+( <TRATC>-1,0)*THETA< MH-I ) 1*44) 

1102  CONTINUE 

00  1103  MM=I4,I5 
SUM=SUM+CAPLAM(M,MMJ 

1103  CONTINUE 

SIJM=SUM+CAPL  AM(M,  !6>*f  TRATC**M 
OUEOOT ( M)=  RN*SUM 
0  IFF  (  M  )=  OUEOOT  I M )-  QL.  AST  ( M ) 

1001  CONTINIJF 

C  PROGRAM  LOGIC  FOR  CHECKING  ON  CONVERGENCE  OF  ITERATIVE  METHOD  FOLLOWS 
1L00K  =0 

151  I LOOK  =  I LOOK  +  1 

IF  (  I  LOOK  -  13)  152,152,1540 

152  IF! DABS  I  0 IFF ( I  LOOK ) )-DABS I OUEOOT ( I  LOOK )*.001 ) ) 1 51 , 1 51 , 1 53 

153  CONTINUE 

OSTAR  =  RN*  I  (EPSKl.  >+EPS2(l.  >)/2.  )*  (  I  ( THETA  IK  )*  ITRATC-1.0) 

1+1. >**4)-(TRATC)**4) 

150  CONTINUE 

O********************************************** 

C  NUMERICAL  INTEGRATION  TO  EVALUATE  2H2I1) 

CONST  =  (1./CII/3. 

Y1N= ( I OUEOOT ( I )-( OUEOOT (1-1) -OUEOOT I  I ) 1/2. )*( ( 1.+DY10X ( 1 . 1**2 )** 
1.5 >  +  ( OUEOOT I  1 3 >  +  (OUEOOT ( 1 3)-QUED0T 1 13-1 ) 1/2. )*( ( l.  +  (0Y20X( 1. 1**21 
2 ) ** . 5 ) ) * ( 1 .-GMS ) 

EVEN=0.0 
0l)D=0.0 
ISTOP  =  1-2 
NC=N 

00  120  KK  =2, ISTOP, 2 
CKK=KK 

EVEN  =  E.VEN+  ( ( ( OUEOOT ( KK)+QUED0T ( KK+1 1 1/2. J*( II, ♦  DY10X(CKK/CI >** 
12  >  ** .5 )♦( ( OUEOOT (KK+I >+QUEDQT( KK+ 1*1 )  1/2  .  )  *  ((1.+  DY2DX(CKK/C I )** 
2?  )  **.5  )  >  *  ( (  (CKK/CI  )**3.  )  -GMS*  (  (CKK/CI  1**2  )  ) 

CKK0=CKK+1 . 

000  =  000  +  (((  OUEOOT  ( KK+1 )+  OUEUOTl KK+2 11/2. )*( ( 1.+0Y10XI CKKQ/C 
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1 1 >**2  >**.S>  +  {  (CHIFOOTIKK+l  +  I  >+C)UFD0TI«  +  I+?n/2.  »*I  1 1  .  +  L1Y20X I CKKI1/C 
?!  >**2  >**.•>>  >*<  ( ICKKU/CI  >**3I-GmS*<  ICKKII/CI  )**?  ) ) 

120  CONTINUE 

Y1L={  I  lOIJFOOT  1 1  l+QUEUOT  (?>)/?,  >*(  ( 1.+UY1DXI  1. /Cl  )**2 ) **. 5  I +  (  < 

imjEortT  n+i  i+oufitot  n+? » >/2. >*i  i  i.+oy2dxi  ucn**?)**.5ii« 

?(  u./cn**3.o-  gms*i i i./ci >**2> > 

C  APPLICATION  DP  SIMPSONS  RULE  GIVES  FOR  ZH2I1)  THE  FOLLOWING 
ZH2(l)  =  CONST  *(4.*(Yll+  00D»+2.*EVEN  +  YIN) 

C  NUMERICAL  INTEGRATION  TO  EVALUATE  ZH2IK)  K=  2.N 

00  1?1  KEXP=2,NC 
C  YK(ZERO)ANi)  Y<  I  N )  =  ZERO 

000=  0.0 

EVEN  =  0.0 

00  122  KK=  2.IST0P.2 
CKK=KK 
0=CKK/C I 
CKK0=CK<+1. 

On=CKKO/CI 

AKEXP  =  KEXP 

EVEN  =  EVEN  +1(1 OUEUOT I  KK )  +  OIJEOOT (KK+1 ) 1/2. )*  | (1.+  DY1DXI0) 
l**2)**.5)+( ( OUEOOT ( KK+ I )+  QUEDOTIKK+l+I )>/2. )* H 1 . +DY2DX I  0 ) 
2**2)**.5))*lin«l  AKEXP+l.>-2.*ID**(AKEXPn  +  D**(AKEXP-l.  I  )  J 
000  =  000  +1(1 OUEOOT I KK+1 )+  OUFOOT I KK+2 )) /2 . ) *  I  1 1 ,  +  OYlOX I  Of)  1**2 ) 
1**.5>+  (I  OIJEDOT  IKK+l+I)  +  OllFOOTIKK  +1  +  2)  1/2.  1*1  1 1 . +0Y2DXI DO  1**2 
2>**.5>>*<00**IAKEXP+1.>-2.*(0D**AKEXP)+D0**(4KEXP-1. I J 

122  CONTINUE 

YK 1= I  I  I OUEOOT 1 1 )+QUF|10T (2))/2.)*( 1 1 ,+OYlOX ( 1 , /C I)**2  >**. 5 1  + 

II  IOUEOOTI  I  +  D+OUEDOTI 1+2  >  )/2.  )  *  I  1 1.+0Y2DX  1 1. /C I )**2 )** . 5 1 ) * 

2  t I l./CI >**IAKEXP+1. )-2.*I I l./CI >** AKEXP )  +  I 1 . /C  I  )** I AKEXP-l . J ) 

C  APPLICATION  OF  SIMPSONS  RULE  GIVES  FOR  ZH2IK)  THE  FOLLOWING 

ZH2!KFXP)=  CONST*  1 4.* 1 YKl+ODO >  +  2 ,*EVEN I 

121  CONTINUE 

C  EVALUATION  OF  THE  COLUMN  MATRIX  ELEMENTS  FOLLOWS 
00  123  KK=1,NC 

CnLUMN(KK>=  ZH(KK>  -I0STAR*ZH1(KK))/AKX  +  ZH2IKK) 

123  CONTINUE 

C***» ********************* «***««** t********************************** 

C*  a**************  *********** t******************* 

EVALUATION  OF  THE  A(H)'S  FOLLOWS 
KGO  =  K  ,  IGO  =  I 
00  127  KGO  =1 *N 
AIKGO 1=0.0 
00  127  IGO  =1,N 
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A(KGO)=A(KGO>+COLUMN<  IGO)*Z(KGO,  IGO) 

127  CONTINUE 

C** **************************** ************************************** 

c******************************************************************** 

C  EVALUATION  OF  CAPPSI(MO)  MO=l,K 
CONST?  =8S*TGS  -  QST AR/AKX 
00  129  M0=1.K 
CKK=K 
AMO=MO 

EX (MO )=  (?.*AM0-l.l/(2.*CK!O 

CAPSI1  =  A(1)*(EX(M0)-GMS)MEX(M0>**2. ) 

CAPPS  I ( MO  >  =  CAPS  1 1 
00  12A  MOO  =2,N 
A  MOO  =  MOO 

CAPPS  II  MO >=  CAPPS  I (MO >♦  A(M00)*(EX(M0I**( AMOO+l. >-2. *( EX ( MO) **4M00 
1)4-  EX(MO)**f  AM00-1.>) 

128  CONTINUE 

C******************************************************************** 

C  EVALUATION  OF  THETA  (I) 

THETA (MO)*  CAPPSI ( MO )+  C0NST2  *  ( EX ( MO )**2-EX( MO) ) 
TBAT(MO)=THETA(MO)*(TRATC-1.0)+1.0 

129  CONTINUE 

C****************************************************** 

£********«************************************************ 

C  CHECKING  FOR  DIVERGENCE  OF  SOLUTION 
00  12R2  MA=1,K 

IF  (THETA(MA)-H>1280, 1280, 1281 

1281  IF  (THETAIMAI-2, 0)1282, 1282, 1283 
1280  THETA(MA)=-1.0 

GO  TO  1282 
1283  THETA (MA  )  =  2.0 

1282  CONTINUE 

C********** *********************************************** 

C********* *************************** ********************* 

C  AVERAGING  THETA  FOR  THF  NEXT  ITERATION 
P COUNT® ICOUNT+1 
00  1285  NIIB=1 ,  I 

THET  A(NUB)®(T6MPAVINUB)*( FCOUNT— 1 .0>+THETA( NUR ) ) /FCOUNT 
TFMPA V ( NUfl  >  =THETA ( NUB ) 

1285  CONTINUF 

C  ******  ***  **  *  *  *  *  *  *  **  * *****  **  *******  **  *  *  *  **  **  ****  *  *** *  * 

C*********************** ****************** ************* 

C  THE  NFXT  COMPUTATION  SAVES  THE  N-TH  OUEOOT  ITERATION  BEFORE  GOING  ON 


82 


A  EDC-TR-69-88 


no  l ao  in  =1,13 

OLAST(IO>=  OUEOOT ( 10) 

140  CONTINUE 

(;*•***«****************«***«*******»*********#*#«**«**« 

C  ICOUNT  IS  THE  COUNTER  FOR  KEEPING  TRACK  OF  THE  NUMBER  OF  ITERATIONS  PERFORMED 
ICOUNT  s  ICOUNT  1 
IF  I  I COUNT -30 >2000 *2 000, 154 

C*4**»**i|c*S**$«4**«**********4******«<t«****:er*********,|i*«*********«***:t 

1540  WRITE(6,603)  ICOUNT 

603  FORMAT  ( • 1THE  ITERATION  CONVERGED  ON  THE  t'I2'»TH  ITERATION.  RESOl 
ITS  FOLLOW. • ) 

C******«******************«*****«**************4* ***************** *44* 

C  CALCULATION  OF  THE  FIN  PERFORMANCE  PARAMETERS  FOLLOWS. 

C  QRSUM  REPRESENTS  THE  ACTUAL  FIN  HEAT  TRANSFER  BY  RADIATION. 

154  QRSUM  =  0.0 

00  190  LA  =1,1 
FLA=LA 

QRSUM  =  QRSUM  C  OUEDOT I L  A  > /Cl  >  *  I  <  1  .♦(  DY10X  (  I  FL  A-.  5  >/C  I  >1**2  )**.  5  > 

QRSUM  =  QRSUM  + 1 QUEOOT ( L A+ 1 ) /C I >* HI .  + 1 DY2DX I  I FLA-.5 ) /C I>)**2 >** 

1.5) 

190  CONTINUE 

QRSUM  =  QRSUM  +  DSTAR  *  ( Y1XI 1. >-Y2X( 1. )) 

C  OCSUM  REPRESENTS  THE  ACTUAL  FIN  HEAT  TRANSFER  BY-  CONVECTION. 

QCSUM  =  0.0 
00  191  LB  =  1 , 1 
FLB=LB 

QCSUM=  QCSUM  +  (B/CI>*ITHETA(LB)-l.>*(HlXHFLB-.5)/CI )*( I 1.+(DY1DX( 
llFLB-.5)/Cm**2l**.5»+  H2X( <FLB-.5)/CI  )  *  (( 1 .  + 1 DY20X  H  FLB-.  5  )  /C I )  I 
2**2)**. 5) > 

191  CONTINUE 

OCSUM=  QCSUM  ♦  B*  <  THETA ( I >“l.)*C (H1X ( 1 .  )+H2X{ 1.))/2.)*(Y1X(1.)- 
IY2XI1.)  1 

C  ORBASE  REPRESENTS  THE  ACTUAL  BASE  HEAT  TRANSFER  BY  RADIATION 
ORBASE  =0.0 

ARE6R  =  BASE  - ( Y1 X 1 0.0 ) -Y2X ( 0 .0 ) ) 

DO  192  LC  =  I4, 1 5 
F IB=  1 B 

ORBASE  =  OR0ASE  ♦  QUEOOT I LC )*  (AREAR/FIRI 

192  CONTINUE 

C  OCBASE  REPRESENTS  THE  ACTUAL  BASE  HEAT  TRANSFER  RY  CONVECTION. 

OCBASE  =-B*( (H1X(0.0)+H2X(0.0> )/2. )*AREAB 
C  ORB  10  IS  THE  IOFAL  BASE  RADIATION  HEAT  TRANSFER. 

0RB10=RN*(1.-TRATC**4)*AREAB 
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c  ocnro  is  the  iofal  base  convection  heat  transfer 
OCBin  =  OCBASE 

C  ORFID  IS  the  iofal  fin  radiation  heat  transfer. 
ciRPin  -o.o 
no  193  LD-1 *  I 
FLD-LO 

ORF If)=OR FID+RN*(1,-TAATC**4)*| ( 1 ,+ ( DY10X ( <  FLO-. 5 ) /C 1 >**?)**. 5  l/C !♦ 
1  (1.+  IOY2IIXI  IFLO-.S>/CI >**2>**.5>/CI > 

193  CONTINUE 

nRFID=0RFI0+AN*<l.-TRATC**4>*( Y1X( l.»-Y2X(  I.  ) ) 

C  OCFID  IS  THE  IDEAL  FIN  CONVECTIVE  HEAT  TRANSFER. 

0CFID=0.0 
DO  19*  LF  =  1,1 
fle  =LF 

OCFID  =OCFin-B*(HlX((FLE-.5>/CI >*(  ( 1.  +  DY10X ( I FLE-. 5  I /C I  )**2>**.5>  + 
1H2X< (FLE-.5)/CI>*( (1.+DY2DXI (FLE-.51/CI )**2I**.5)  1*11. /Cl ) 

194  CONTINUE 

OCFin  =0CFID-8*((H1X<1.)+H2X(1.)>/2.)*(Y1X(1.I-Y2X(  l.»> 

ETA=  (  ORStJM  +GCSUM)/  (ORFIO  +  OCFID) 

APEM I T=  ( ORSIIM  +ORBASE  -OSTAR*  I Y1X  ( 1.0I-Y2XI 1.0) ) )/(RN*( 1.4 
1TRATC**4  >*3ASE ) 

FFLUX=(QRSUN+ORBAS6+QCSUM+QCBASE>/l I0RBID+0CBI0)*(BASE/AREA8) ) 
WRITE  (6.641)  ORSUM  ,OCSUM, ORBASE, OCBASE. ORFI 0, OCFID,  ORMO.OCBID 

641  FORMAT  (//•  ACTUAL  FIN  RADIATIVE  HEAT  FLUX  =  'E13.5.10X,1  ACTUAL  F 
1IN  CONVECTIVE  HEAT  FLUX  a  'E13.5,/,*  ACTUAL  BASE  RADIATIVE  HEAT  FL 
?U*  =  ,F13.9»9X.i  ACTUAL  BASE  CONVECTIVE  HEAT  FLUX  =*S14.5,/,*  IDEA 
3L  FIN  RADIATIVE  HEAT  FLUX  =  'E13.5.I1X.'  IDEAL  FIN  CONVECTIVE  HF AT 
4  FLUX  =  *613.5,/,*  IDEAL’ BASE  RADIATIVE  HEAT  FLUX  =  'E13.5.10X, 

5*  IDEAL  BASE  CONVECTIVE  HEATFLUX  *  'E13.5) 

WRITE!  6.640  )'  ETA.APEHIT  .EFLUX 

640  FORMAT ( *  FIN  EFFICIENCY  *  «E13.5,/t*  APPARENT  EMITTANCE  -  *613.5 
1./,*  EFFECTIVE  HEAT  FLUX  =  *E13.5> 

155  DO  13P  IT-1, 1 
NT=1 T  +1 

WRITE  (6,604)  IT  ,OUEDOT( I T > , I T.DI FF( IT) , NT,OUEDOT( NT  1 , NT, DI FF ( NT ) 

604  FORMAT  <•  OUEDOT (  *  12 • >  =  • E14. 7 ,4X, *DI FF( ■ 1 2 • )  -  IE14.7******  OUE 
1D0T ( *  1 2 1 )  =  *E14.7,4X,*DIFF(* 12* )=  *E14.7) 

138  CONTINUE 

WRITE! 6,48 ) (OUEDOT (NU) *NUaI4, 16 ) 

48  FORMAT  1 6E20 .5 ) 

DO  139  I  A  =  1 , N 
WRITE(6,605)IA,A(  IA) 

605  FORMAT ( *  A(*!2*)  =  *614.7) 
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139  CONTINIJF 
ihalf  =1 n 
I'U  1^1  IT  =  1,IH4LP 
IS=|T+IH4I> 

WKI I H 6, 62 3) IT, THETA!  IT),  IT,  THAT  (  IT)  ,  I S , THFTA ( I S  )  , I S , TR AT <  IS) 

623  FORMAT!'  ThFTAI'I?')  =  'E14.7,4X, 'TRAT! ' 12*  )  =  •  E 14.  7  * *****  THFTA  ( 

1>I?>>  =  ‘hlA.y.AX,  •  TRAT ( • I  2 1  )  =  *  F 1 4. 7 ) 

141  CONTINUE 

WR  I  T  r  I  6, 77ft  >ft 

77ft  FORMAT ( •  THE  CONVECTION  NUMRFR= • F6. 2 ) 

WR I T  F ( 6, 777  >  RN 

777  FORMAT ( •  THF  R  AO  I  AT  ION  NUMftER= • F6. 2 > 

WRI1 F(6,779)BASE 

779  FORMAT ( •  THE  FIN  PITCH  TO  LFNGTH  RATIO  IS  'F6.2I 
WRITE  ( 6  *  780 ) EF  I N 

780  FORMAT! ■  FIN  SURFACE  FMISSIVITY  =  'FA. 4) 

49  CONTINUE 

50  CONTINUE 

51  CONTINUE 
GO  TO  135 

156  WRITE  (6,626) 

626  FORMAT  C1THE  ITERATION  SCHEME  FAILED  TO  CONVERGE  IN  30  ATTEMPTS  A 
1NO  HAS  BEEN  TERMINATED.  VALUES  EXISTING  AT  TERMINATION  FOLLOW.') 

GO  TO  154 
135  STOP 

******************************************************************************** 


COMMENTS  On  SUBROUTINES 


SUBROUTINE  VERTDR  IS  AN  AEUC  IRM  360/50  LIBRARY  SUBROUTINE  MODIFIEO  TO  RELAX 
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INVERSION  SUBROUTINE. 

******************************************************************************** 

C 


85 


noooon  *0000000 


AEDC-TR-69-88 
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surroutinf  factor  is  a  double  precision  surrootinf  fur  cai.Cim.at img  the 
SHAPEFACTORS  OF  THE  V/AK IOUS  FI.ENFImTS  USED  1 N  THF  RADIATION  ANALYSIS.  FACTOR 
IS  A"1  AMPLICATION  OF  THE  STRING  POLYGON  METHOD  OF  CALCULATING  SHAPFFACTORS. 
FACTOR  WAS  I'FVFLOPFO  HV  OR.  HALTER  FROST  OF  THE  UNIVERSITY  OF  TENNESSEE 

spacf  institute  ano  hill  soon  re  purlishfo  as  an  afdc  tor. 


SUBROUTINE  VF«TO  IS  AN  AFOC  I RM  360/50  LIBRARY  SUBROUTINE  FOR  OOIIRLF 
PRECISION  MATRIX  INVERSION. 
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